Active learning accelerates ab initio molecular dynamics on reactive energy surfaces

•Autonomous acquisition of datapoints for reactive force field training•Neural nudged-elastic band•Neural molecular dynamics•Switching from entropic to thermodynamic intermediate in solvents with high polarity In silico elucidation reaction mechanisms using density functional theory (DFT) can explain and predict experimental observations. However, because the computational cost accurate DFT simulations, theoretical studies are often restricted systems fewer than 100 atoms at a few stationary points on potential energy surfaces. These be insufficient describing challenging where dynamical effects important. Herein, we report low-cost transferable pipeline that accelerates ab initio dynamics by factor 2,000. It consists high-throughput computation, active learning, transfer learning train high-quality fields based neural networks. reproduce underlying surface enable simulations. We anticipate our will lead affordable mechanistic complex, experimentally relevant systems. Modeling chemical reactions typically requires (AIMD) simulations due breakdown transition state (TST). Reactive AIMD limited lower-accuracy electronic structure methods weak statistics quantum mechanical energies forces must evaluated femtosecond time resolution over many replicas. data-driven allows treatment same level overall as TST approaches. High-throughput calculations autonomous data coupled graph convolutional neural-network interatomic potentials, allowing inexpensive accuracy. demonstrate approach accurately simulating post-TS three distinct pericyclic reactions, including trispericyclic complex bifurcating surface. This is broadly applicable understanding predicting outcomes large, previously intractable Transition (TST) lies heart most chemistry approaches reactivity. Reaction rates product selectivities quantified, thus predicted, free differences between reactants states (TS), whose atomic configurations taken static (PES). Effective combine numerical optimization intuition, such band (NEB) eigenvector following (EV), routinely used determine geometry TSs near-experimental Typically, geometries corrections, high-level correlated CCSD(T) energies. The nature TST, however, fails effects, post-transition-state bifurcation,1Hare S.R. Tantillo D.J. Post-transition bifurcations gain momentum – current field.Pure Appl. Chem. 2017; 89: 679-698Crossref Scopus (91) Google Scholar, 2Tan J.S.J. Hirvonen V. Paton R.S. Dynamic intermediates radical cation Diels-alder cycloaddition: lifetime suprafacial stereoselectivity.Org. Lett. 2018; 20: 2821-2825Crossref PubMed (10) 3Yang Z. Jamieson C.S. Xue X.-S. Garcia-Borràs M. Benton T. Dong X. Liu F. Houk K.N. Mechanisms involving intermediates.Trends 2019; 1: 22-34Abstract Full Text PDF (18) Scholar which responsible distribution transformations. include SpnF catalyzed biosynthetic pathway Spinosyn A4Yang Yang S. Yu P. Li Y. Doubleday C. Park J. Patel A. Jeon B.S. Russell W.K. H.-W. et al.Influence water enzyme energetics ambimodal [6+4]/[4+2] cycloaddition.Proc. Natl. Acad. Sci. USA. 115: E848-E855Crossref (48) several LepI enzyme.5Ohashi Hai Chen Tang M.C. Sato Watanabe K. SAM-dependent enzyme-catalysed natural biosynthesis.Nature. 549: 502-506Crossref (94) study these (AIMD).6Iftimie R. Minary Tuckerman M.E. Ab dynamics: concepts, recent developments, future trends.Proc. 2005; 102: 6654-6659Crossref (130) this approach, nuclei treated Newtonian bodies, problem solved fly (QM) femtosecond-scale step. hindered compounding costs long trajectories (102Tan 4Yang QM evaluations) multiple replicas Scholar,3Yang Scholar). Studies small basis sets. also prevents comprehensive exploration different conditions, substituents. Recent advances machine (ML) techniques modularization network (NN) automatic differentiation routines7Abadi Agarwal Barham Brevdo E. Citro Corrado G.S. Davis Dean Devin al.TensorFlow: large-scale heterogeneous distributed systems.arXiv. 2016; https://arxiv.org/abs/1603.04467Google facilitate interplay physical sciences statistical learning.8Butler K.T. Davies D.W. Cartwright H. Isayev O. Walsh Machine materials science.Nature. 559: 547-555Crossref (1157) Scholar,9Wang W. Axelrod Gómez-Bombarelli Differentiable control learning.arXiv. 2020; provide an opportunity speed up costly Born-Oppenheimer (BOMD) orders magnitude. NN-based potentials built variety representations retain accuracy training executed comparable classical fields. To date, have been developed, starting seminal Behler-Parrinello approach10Behler Parrinello Generalized representation high-dimensional potential-energy surfaces.Phys. Rev. 2007; 98: 146401Crossref (1507) ANI,11Smith J.S. Roitberg A.E. ANI-1: extensible cost.Chem. 8: 3192-3203Crossref (645) 12Smith Nebgen B.T. Zubatyuk Lubbers N. Devereux Barros Tretiak Approaching cluster general-purpose through learning.Nat. Commun. 10: 2903Crossref (185) 13Devereux Smith K.K. Extending applicability ANI deep sulfur halogens.J Chem Theory Comput. 16: 4192-4202Crossref (33) tensor network,14Schütt Arbabzadah Chmiela Müller K.R. Tkatchenko Quantum-chemical insights networks.Nat. 13890Crossref (585) SchNet,15Schütt Sauceda H.E. Kindermans P.J. SchNet - A architecture molecules materials.J. Phys. 148: 241722Crossref (488) Scholar,16Schütt Kessel Gastegger Nicoli K.A. SchNetPack: toolbox atomistic systems.J. Theor. 15: 448-455Crossref (116) net,17Zhang L. Han Wang Car E Deep scalable model mechanics.Phys. 120: 143001Crossref (427) MEGNet,18Chen Ye Zuo Zheng Ong S.P. Graph networks universal framework crystals.Chem. Mater. 31: 3564-3572Crossref (198) DimeNet,19Klicpera Groß Günnemann Directional message passing graphs.arXiv. : 03123Google OrbNet.20Qiao Welborn Anandkumar Manby F.R. Miller T.F. OrbNet: symmetry-adapted atomic-orbital features.J. 153: 124111Crossref (43) effectively learn short-range many-body interactions all-atom coarse-grained potentials.21Wang Coarse-graining auto-encoders dynamics.npj 5: 125Crossref (53) Scholar,22Ruza Schwalbe-Koda D. Harris W.H. Temperature-transferable coarse-graining ionic liquids dual networks.J 164501Crossref (14) their reviews, Noé23Noé Clementi simulation.Annu. 71: 361-390Crossref (161) Jensen24Coley C.W. Eyke N.S. Jensen K.F. Autonomous discovery part I: progress.Angew. Int. Ed. Engl. 59: 22858-22893Crossref (80) covered important aspects accelerating ML techniques. Particularly area catalysis solid-state diffusion, proposed accelerate location TS minimum paths.25Yao Herr J.E. Toth Mckintyre Parkhill TensorMol-0.1 chemistry: augmented long-range physics.Chem. 9: 2261-2269Crossref 26Koistinen O.P. Ásgeirsson Vehtari Jónsson Nudged elastic accelerated Gaussian process regression inverse distances.J. 6738-6751Crossref (23) 27Garrido Torres J.A. Jennings P.C. Hansen M.H. Boes J.R. Bligaard Low-scaling algorithm nudged surrogate model.Phys. 122: 156001Crossref (61) Similarly, simulation aimed properties single organic molecules28Hachmann Olivares-Amaya Atahan-Evrenk Amador-Bedolla Sánchez-Carrera Gold-Parker Vogt Brockway A.M. Aspuru-Guzik harvard clean project: screening design photovoltaics world community grid.J. 2011; 2: 2241-2251Crossref (382) 29Lin Beh E.S. Tong Q. Valle Aziz M.J. Gordon R.G. redox-flow battery alloxazine-based electrolyte.Nat. Energy. 16102Crossref (280) 30Gómez-Bombarelli Aguilera-Iparraguirre Hirzel T.D. Duvenaud D.K. Maclaurin Blood-Forsythe M.A. Chae H.S. Einzinger Ha D.-G. Wu al.Design efficient light-emitting diodes virtual approach.Nat. 1120-1127Crossref (465) 31Pun G.P.P. Batra Ramprasad Mishin Physically informed artificial modeling materials.Nat. 2339Crossref (78) increasingly targeting reactivity.32Wei J.N. Neural prediction reactions.ACS Cent. 725-732Crossref (229) 33Ioannidis E.I. Gani T.Z.H. Kulik H.J. molSimplify: toolkit automating inorganic chemistry.J. 37: 2106-2117Crossref (67) 34Friederich dos Passos Gomes G. De Bin Balcells dihydrogen activation space surrounding Vaska’s complex.Chem. 11: 4584-4601Crossref 35Grambow C.A. Pattanaik Green energies.J. 2992-2997Crossref Here, show combined application (HT) use learning,36Smith B. Less more: sampling learning.J. 241733Crossref (220) net methodology investigate effect solvent viscosity concertedness acquiring sparse data, deploying expansive MD define concerted those gap < 60 fs bond changes,3Yang First, two were studied: classic Claisen rearrangement allyl vinyl ether (Scheme 1A) Diels-Alder acrylaldehyde cis-butadiene 1B). was then tested complicated 8,8-dicyanoheptafulvene 6,6-dimethylfulvene 1C), behavior cannot explained TST.37Liu C.Y. Ding S.T. Cycloadditions electron-deficient 8,8-disubstituted heptafulvenes electron-rich 6,6-disubstituted fulvenes.J. Org. 1992; 57: 4539-4544Crossref Scholar,38Xue X.S. Ambimodal Trispericyclic dynamic periselectivity.J. Am. Soc. 141: 1217-1221Crossref (29) revealed unique each would difficult elucidate traditional cost. Specifically, reaction, conclude relationship bifurcation heavily dependent dielectric constant charge-separated intermediate. envisage produce system-specific faithful trained high-accuracy data. empower intractable. key components loops (1) HT rapid, parallel data; (2) networks; (3) rapid inference sample configurational versions NEB 2). Because uncertainty not necessary AL it categorized membership query synthesis, Settles.39Settles Active literature survey.http://digital.library.wisc.edu/1793/60660Date: 2009Google For constructed independent cycles acquire resulting reaction-specific potentials. cycle, 14,976 7,927 M062X/def2-SVP performed, respectively. more steps leading product. converged total 31,397 obtained theory. All attained sub-kcal mol−1 corresponding test sets (Table S2). further ascertain models, relative quantities 1). mean absolute error (MAE) 0.4 kcal mol−1, attaining sub-chemical reference M06-2X/def2-SVP method. Despite explicitly Hessians, calculate various quantities, MAE Gibbs correction term, Gcorr, had reasonable value 2.3 mol−1. majority deviation came Scorr, known large errors presence low-frequency vibrations. Several quasi-harmonic corrections could applied better estimates.40Ribeiro R.F. Marenich A.V. Cramer C.J. Truhlar D.G. Use solution-phase vibrational frequencies continuum models solvation.J. 14556-14562Crossref (590) Scholar,41Grimme Supramolecular binding thermodynamics dispersion-corrected theory.Chemistry. 2012; 18: 9955-9964Crossref (894) ScholarTable 1Relative NN [in kcal/mol]ClaisenSpeciesRel. EaM06-2X/def2-SVP respective reactants.ZPEHcorrTScorrGcorrbGcorr = ZPE + Hcorr -TScorrCL-R0.0 (0.0)73.7 (74.6)5.0 (5.2)25.5 (24.2)53.2 (55.4)CL-TS32.9 (32.4)71.0 (74.0)4.2 (4.8)24.0 (22.3)51.2 (55.9)CL-P−16.6 (−16.6)72.9 (74.5)5.0 (5.1)24.7 (24.3)53.2 (55.2)Diels-AlderDA-R0.0 (0.0)95.5 (93.4)6.0 (7.3)27.2 (29.6)74.3 (71.1)DA-TS15.2 (14.3)94.4 (95.0)6.4 (5.7)28.4 (25.9)72.4 (74.8)DA-P−49.1 (−49.2)96.2 (98.2)5.4 (5.2)25.3 (24.9)76.3 (78.5)TrispericyclicTP-R0.0 (0.0)180.8 (181.1)12.7 (12.5)56.5 (55.0)136.9 (138.7)TP-TS18.4 (7.6)181.1 (181.9)12.2 (11.9)43.9 (39.1)149.6 (154.7)TP-P1−19.5 (−19.6)184.4 (184.4)11.6 (11.6)39.5 (38.5)156.5 (157.4)TP-TS2a5.6 (5.4)182.3 (182.6)11.5 (11.4)39.3 (37.9)154.4 (156.0)TP-P2a−27.2 (−27.4)184.3 (185.1)11.3 (11.1)38.7 (37.2)156.9 (159.0)TP-TS2b6.1 (6.2)182.9 (182.6)11.3 (11.5)38.0 (38.3)156.2 (155.8)TP-P2b−23.4 (−24.1)184.4 (183.7)11.8 (11.7)40.9 (39.0)155.2 (156.5)ZPE, Hcorr, Gcorr abbreviations zero-point energy, enthalpy (without ZPE), entropy, Species defined Scheme 1. T 298.15 Ka reactants.b -TScorr Open table new tab ZPE, K efficiently upgrade solvation leveraged (TL) strategy36Smith larger set (def2-TZVPD) SMD (for details see Supplemental information). common (cyclohexane, chloroform, acetone) about 10%–20% additional S1). Furthermore, TL performed double-hybrid vacuum, DLPNO-DSD-PBEP86-D3BJ/def2-TZVPD (abbreviated DH/def2-TZVPD), role functional. transfer-learned out-of-sample (Figures S1–S5 Claisen, Figures S6–S10 Diels-Alder, S11–S15 trispericyclic). initiated vicinities QM-optimized levels theories forward backward directions (see Computational details). sections, mainly discuss results (SMD-)M06-2X/def2-TZVPD order focus dynamics. foundational higher quality DH/def2-TZVPD vacuum comparisons. utilize Langevin thermostat friction coefficient total, 500 pairs vicinity CL-TS solvent. Of these, number valid (those resulted CL-P one end CL-R other) ranged 384 474 S3). calculated breaking C–O formation C–C γ,δ-unsaturated aldehyde quantitative measure reaction. breaks (C–C forms) when length exceeds (falls below) 1.7 Å. value, while somewhat arbitrary, consistent convention adopted previous studies.3Yang histograms gaps (breaking) (forming) lengths (SMD-)M06-2X/def2-TZVPD-optimized presented Figure Regardless environment, percentage dynamically (with less fs) 100% cyclohexane, 97% 98% chloroform acetone, Both median fall narrow range 26–36 fs. note distributions acetone spread out increased stability transient species (after but before forming) increase dynamics(NRMD) showed no significant terms (97%), (37 35 fs, respectively). list relating Table S3. NRMD 400 generated suggest there subtle concertedness. polar 90% concerted, 18, 20, 14 picoseconds (ps), although slightly bimodal observed another peak around 200 ps (Figure case decreased 92% 57% S3) much broader clearly 86 ps. disentangle polarity, all four media fixed values 10−4, 10−3, 10−2 a.u. trends across arose S4). reasons longer magnitude imaginary frequency DA-TS SMD-M06-2X/def2-TZVPD −468.57, −484.73 −481.50 cm−1 respectively, suggesting flatter cyclohexane possible observed. Comparing NRMDs M06-2X/def2-TZVPD 93% 89%, 28 43 18 31 possibly inclusion MP2 correlation explicit Grimme’s D3 dispersion term,42Grimme Ehrlich Goerigk Effect damping function corrected theory.J. 32: 1456-1465Crossref (9713) Scholar,43Grimme Antony Krieg parametrization (DFT-D) 94 elements H-Pu.J. 2010; 132: 154104Crossref (22221) stabilize after Experimental room temperature yield 0: 1.17: 1 TP-P1[4+6]: TP-P2a[6+4](+tautomer): TP-P2b[8+2] day. Recently, al.38Xue carried 500-fs BOMD ωB97X-D/6-31G(d) TP-TS1, 87% 142 TP-P1, 3% TP-P2a TP-P2b. TP-TS2a TP-TS2b, reported result TP-P2b, From this, concluded products TP-P2b experiment control, connection made yield. understand choice method system, simulated medium significantly timescale 5 direction, so account lifetimes both friction. addition connecting minima, “recrossing” “hovering” “recros